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I. Summary 


Analysis of interplanetary trajectories is a crucial area for both manned and 
unmanned missions of the Space Exploration Initiative. A deep space maneuver 
(DSM) can improve a trajectory in much the same way as a planetary swingby. 
However, instead of using a gravitational field to alter the trajectory, the on-board 
propulsion system of the spacecraft is used when the vehicle is not near a planet. 
There are occasions (broken plane maneuvers are one example) where the 
advantages gained at the endpoints of the trajectory outweigh the cost of the DSM. 
The purpose of this study is to develop an algorithm to determine where and when 
to use deep space maneuvers to reduce the cost of a trajectory. The approach taken 
to solve this problem uses primer vector theory in combination with a non-linear 
optimizing program. Primer vector theory applies the calculus of variations to the 
trajectory problem in order to minimize AV. A set of necessary conditions on a 
Lagrange multiplier called the primer vector arises from the analysis, and this 
primer vector indicates whether a deep space maneuver will be beneficial. Deep 
space maneuvers are applied to a round trip mission to Mars to determine their 
effect on the launch opportunities. Other studies which were performed include 
cycler trajectories and Mars mission abort scenarios. It was found that the software 
developed was able to quickly locate DSM's which lower the total AV on these 
trajectories. 
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Abbreviations 
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see equation A. 10 
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index ranging from 1 to m, or from 1 to 3, as specified 

index ranging from 1 to p 

number of controls 

number of states 

number of constraints 

number of final states which are fixed 
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1. Introduction 


Background 

With the continuing efforts to develop Space Station Freedom, efforts are also 
being directed toward the next steps in our space program; a return to the moon, and 
then a manned expedition to Mars. In relation to these efforts, studies are being 
performed to determine optimal trajectories for these missions. Trajectories can be 
simulated and numerically optimized on a computer. These simulations can 
determine the encounter dates which minimize some cost function, say initial mass 
of the vehicle or total AV. However, other missions with different architectures and 
lower costs may exist. The addition of a planetary swingby or a deep space 
maneuver (DSM) changes the architecture of a mission and may reduce the overall 
cost of the mission. 

Since the location of a swingby is determined by the position of the swingby 
planet, time is the only independent variable for the addition of a swingby. 

Therefore, to determine if a swingby will reduce the cost of a trajectory, a second 
simulation can be performed with the swingby architecture. By varying the time of 
the swingby, the minimum cost of the new architecture can be determined and 
compared to the cost of the original trajectory. The criteria for determining if a DSM 
will reduce the cost of a mission is not so simple. One example of when a DSM will 
be beneficial is in the transfer of a satellite between two non-coplanar orbits. Say 
there is a satellite at point A in orbit 1 in Figure 1.1, and we want to rendezvous 
with a second satellite at point B in orbit 2. This rendezvous can always be 
accomplished with a two burn transfer trajectory. The plane of this trajectory is 
determined by the points A and B, and the center of attraction of the orbits. As in 
the figure, this plane may be highly inclined. The out of plane component of the 


velocity on the transfer orbit will require large AV's at points A and B. A DSM can 
be used to reduce these out of plane components, and thus the total AV. If the 
satellite remains in the plane of orbit 1 after the initial burn, an additional burn can 
be applied when the satellite reaches the line of nodes. This burn would be used 
only to change the plane of the transfer orbit from that of orbit 1 to that of orbit 2. A 
final burn can then be applied at Point B to match the velocity of the second satellite. 
Previous studies 1 have shown that allowing some of the out of plane AV's to take 
place at Points A and B can reduce the total AV required even further. The DSM 
would not take place right on the line of nodes, but rather some distance from it. 

The location and time of a DSM are independent, and therefore the one- 
dimensional search used for the swingby is no longer applicable. Since a DSM can 
be located anywhere in space, a four-dimensional search of time and the three 
coordinates of the DSM requires a great amount of computation. However, if a 
reliable initial estimate can be generated for the time and location of a DSM, this 
search can be reduced to a reasonable effort. 

Primer vector theory, introduced by Lawden^, can be used to make these 
initial estimates. Previous work on this subject stemming from the work of Lawden 
has been done by Jezewski and RozendaaP, GlandorH, and Lion and Handelsman^, 
among others. MULIMP^ (Multi-Impulse Trajectory and Mass Optimization 
Program) is a piece of software which was developed at ITT Research Institute, and 
later at Science Applications, Inc. for the Jet Propulsion Laboratory in Pasadena, 
California. MULIMP applies primer vector theory to place DSM's on interplanetary 
trajectories, and may be the tool most commonly used to do so. 

The purpose of this study is to develop new software which uses primer 
vector theory to locate DSM's quickly on interplanetary trajectories. The software 
has been integrated into a trajectory analysis code called IPREP (Interplanetary Pre- 
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Processor) 7 . IPREP previously had the capability to place DSM's on trajectories, 
however the routine which it used to find DSM's was very slow and not reliable. 
This paper provides an introduction to primer vector theory, a method to 
determine the initial estimates for the time and location of a DSM if one will be 
beneficial, and an algorithm to incorporate DSM's on trajectories with swingby's. 
Incorporating a DSM on a trajectory with a swingby presents a difficulty because if 
the DSM precedes the swingby, then the inbound Voo vector to the swingby is not 
known a priori, and the swingby parameters cannot be calculated. The situation is 
similar if the DSM follows the swingby or if more than one DSM is used on the 
same leg of a trajectory. 

Thesis Organization 

The first major section of this thesis is devoted to reviewing all of the 
assumptions and approximations used in the analyses which follow. A description 
of the methods used in IPREP to calculate and optimize trajectories is presented. 

The derivation and assumptions of primer vector theory will then be summarized. 
Finally, the process of incorporating primer vector theory into IPREP will be 
discussed. 

Results are discussed in the next section. In order to verify that the new 
software was functioning properly, a previous study which used DSM's was 
analyzed with the new software. This study was of a Cycler trajectory between Earth 
and Mars. A second study was done on round trip Mars missions to show the affect 
DSM's can have on expanding launch windows. Abort scenarios for Mars missions 
were also examined to determine how DSM's might be used to open more abort 
options for manned missions. Finally, conclusions are drawn and 
recommendations for future work are prescribed. 
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2. Analysis 


Description of IPREP 

IPREP was originally intended for use as an aid in using IPOST 7 
(Interplanetary Program to Optimize Simulated Trajectories). IPOST will optimize a 
trajectory when initial estimates for some of the parameters of that trajectory are 
input. However, the input parameters must describe a reasonable trajectory to 
guarantee convergence in IPOST. IPREP determines reliable initial estimates of the 
Voo vectors and encounter dates to be used for inputs into IPOST. However, IPREP 
can also be used as a stand alone tool for preliminary analysis of interplanetary 
trajectories. 

Inputs to IPREP include the order of the planets to be encountered, 
maneuvers to be performed at each encounter, and a time of flight window from 
each encounter to the next. For each of these windows, there is also an input step 
size which determines how many dates will be checked in that window. For a given 
set of encounter dates, IPREP calculates the trajectory using a patched-conic 
technique. In this technique, the planets are assumed to be point masses, and the 
position and velocity of each planet is found from one of the ephemerides available 
in IPREP 8 . A transfer orbit for each leg of the trajectory is found by solving 
Lambert's problem. With these transfer orbits the Voo vectors are determined. For 
launch or orbit insertion maneuvers, a parking orbit must be defined by the user. 
With the parking orbit defined, the periapsis velocity in the parking orbit is known. 
The energy equation: 


E = 



2 



2 


r P 


( 2 . 1 ) 
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is used to determine I Vperi I / the magnitude of the velocity at periapsis along the 
inbound or outbound hyperbola for orbit insertion or launch, respectively. The 
difference between this value and the magnitude of the velocity in the parking orbit 
at periapsis is then the required AV, with the assumption being that all AV's are 
made to be tangential to the orbit at periapsis. Notice that this method produces a 
magnitude for the AV, but the direction of the AV vector is not determined. 

If the inclination of the parking orbit is defined by the user and is such that a 
plane change must be performed in order to get from the Voo vector to the parking 
orbit, this plane change burn is calculated separately and its magnitude (not a vector 
addition) is added to the AV required for the maneuver. 

Two methods are available in IPREP to determine the AV required for a 
swingby maneuver. The first and more simple method compares the magnitude of 
the inbound and outbound Voo vectors and determines the angle between them by 
taking their scalar product. The angle between the Voo vectors determines the 
required turn angle, 0, for the swingby, which is also given 9 by: 


. 0 1 

sin - = r . 

2 i , r P V ".in 

M- 


(2.2) 


If the required value of rp, the swingby radius, falls between the user input values 
for the minimum and maximum values, then the required swingby AV is the 
difference in magnitude between the inbound and outbound Voo vectors. However, 
if the required rp falls outside of the allowed range, the additional turn angle must 
be made up by the AV. In this case, the AV becomes: 


— y Voo jr, + Voo ^j +|Voo,in V^outl*- 08 ^ ®max) 


(2.3) 
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This method applies the AV after the swingby has been completed, instead of at 
periapsis. 

The second method to calculate swingby AV's uses a non-linear optimizer to 
find the smallest AV which can complete the swingby 9 . The AV is applied at 
periapsis along the inbound hyperbolic trajectory (which in general is not the same 
as periapsis along the outbound hyperbolic trajectory). The non-linear optimizer is 
invoked with the radius of periapsis and the direction of the swingby AV vector as 
controls, and the magnitude of the swingby AV as a cost function. A lower limit of 
1.1 planetary radii is imposed on the periapsis radius to avoid atmospheric 
interaction. 

After all of the AV's have been found, a cost function for that trajectory is 
evaluated 4 : 


cost = 


I 


n=l 


l V ~.in| n 

tol(l,n) 


+ 


V °°,OUt| n 

tol(2,n) 


t AV n | AV DSM>n | tof n 
tol(5,n) tol(6,n) tol(7,n) 


initial mass 
tol(3, 1) 


final mass 
tol(4,N) 


(2.4) 


where N is the number of events in the trajectory, tof n is the time of flight for the 
n lh leg, and tol(l-7,N) is a user defined array of weighting parameters. By weighting 
different parameters in this equation, the user can define the cost function to be the 
initial mass, total AV, or some other function. 

IPREP begins by running the trajectory with all of the dates set to the 
beginning of each encounter window. The last date is incremented and new 
trajectories are run until the end of the last encounter window is reached. At that 
point, the second to last date is incremented, the last date is reset to the beginning of 
its window, and the process continues. The cost function is evaluated for each 
trajectory, and the set of dates which produces the lowest cost function is saved. 
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After all of the combinations of dates has been run, the minimum cost trajectory is 
written to an output file. 

This grid search' method of optimization has advantages and disadvantages. 
The disadvantage is that for trajectories with large encounter windows or many 
events, the number of trajectories to be calculated grows very quickly. For a 
trajectory with N events and m dates to be checked for each event, the total number 
of trajectories to be run is m^. The advantage to the grid search method is that it 
will not become trapped in a local minimum as gradient methods do. IPREP does 
also allow the user to define a maximum total trip time, which may eliminate some 
of the combinations of dates before those trajectories are run. 

Primer Vector Theory 

The derivation of four necessary conditions for optimality on impulsive 
trajectories where the total AV is minimized was presented by Lawden 2 . That 

derivation is summarized in the Appendix of this report. The conditions derived 
are: 

1) The primer vector and its first time derivative are continuous. 

2) During any impulse, the thrust vector must be aligned with the primer. 

3) The magnitude of the primer is a maximum and has a value of one during 
any impulse. 

4) The derivative of the magnitude of the primer is zero at a deep space 
maneuver. 

In arriving at these conditions, the spacecraft was treated as a point mass in an 
inverse square law gravitational field. The thrusts were modeled as impulses, and 
assumed to occur over infinitesimal time duration. Also, the cost function was 
taken to be the total AV, which is given by the rocket equation as: 
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(2.5) 
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M. 


AV = golspl In 

It should be noted that the cost function here is different than the cost function used 
in the grid search optimization in IPREP. However, minimizing total AV is 
equivalent to minimizing the initial mass of the vehicle if no discontinuities in the 
mass of the vehicle are incurred. Mass discontinuities occur if part of the vehicle 
(say a Mars excursion vehicle) is jettisoned. In general, this only occurs at the end of 
a leg of the trajectory, and since primer vector theory will be applied to each leg of 
the trajectory separately, total AV is a good choice for the cost function here. 

It is noted from equations A.20 and A.21 in the appendix of this report, that 
the primer vector, X, satisfies the same differential equation as the 'deviation 
vector' in equation 9.31 of Reference 10. The solution to this equation is given by 

equation 9.45 in Reference 10 as 

= v(t/to) 


Or 


Ov 


-it 


Or 


Ov 


( 2 . 6 ) 


*o 


Therefore, the primer vector can be found at any point along a trajectory from the 
equation: 


T 

= Y(t,t 0 ) 

X() 


t 

*0. 


(2.7) 


where \|/(t,t 0 ) is the overall state transition matrix. If \| / and the initial conditions 
on the primer vector and its derivative are known, then the primer and its 
derivative can be found along the entire trajectory. Assuming we have a mission 
which begins with a launch and ends with an orbit insertion, then the primer vector 
at the initial and final times is known to be a unit vector in the direction of the 
thrust at those points. The parameters of the parking orbits at each end of the 
trajectory provide the information to determine the position and velocity of the 
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spacecraft at periapsis. From this information, and from a Lambert solution 
between the initial and final position of the spacecraft, the direction and magnitude 
of the initial and final AV's can be determined. Equation 2.7 becomes: 


T 

= V(tf/to) 

*0 



>11 

¥12* 

y 

x_ 

t=tf 

>0. 


_V21 

¥22. 



(2.8) 


where the 6X6 matrix \|/(tf,t()) has been replaced with the equivalent four 3X3 
matrices. These leads to 5 : 


X f = + Vi 2^0 


(2.9) 


X { ) = \| f]2 


¥11^0] 


( 2 . 10 ) 


If Y 12 is not singular, we can find the initial conditions on the derivative of the 
primer vector. With these conditions known, the primer vector can be propagated 
along the entire trajectory. 

A problem occurs when we cross into or out of a sphere of influence (SOI). 
The sphere of influence is an imaginary boundary, where inside the sphere all 
motion is assumed to be governed solely by the gravitational attraction of the 
planet, and outside of the sphere, motion is assumed to be governed by the 
gravitational attraction of the sun. To further develop the idea of the SOI, consider 
a spacecraft in space near a planet. Its motion can be modeled as a two body 
problem, where the spacecraft and the planet are orbiting around their common 
center of gravity. The attraction of the sun, which is considerably farther away, can 
be approximated as a disturbing force in that two body problem. As the spacecraft 
moves away from the planet, the ratio of the disturbing force from the sun to the 
force from the planet increases. As we move farther from the planet, the motion 
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can be modeled as the spacecraft orbiting the sun, and the attraction of the planet 
could then be considered the disturbing force. At some distance from the planet, the 
ratio of the disturbing force to the central force will be the same for each of these two 
models. The mean distance at which this is true is defined as the radius of the 
sphere of influence for the planet. The magnitude of this radius will be dependent 
on the direction the spacecraft has moved from the planet. However, if the distance 
between the sun and the planet is much greater than this radius, then the radius is 
approximately constant for all directions^, and the SOI is approximated to be a true 
sphere. The value for the radius of the SOI for each planet is provided by the 
ephemeris information in IPREP 8 . 

Since the initial and final conditions on the primer are determined by the AV 
vectors inside the SOI, the primer vector must first be propagated from its initial 
condition out to the boundary of the SOI. When the primer exits the SOI, Glandorf 
has shown 4 that the primer is continuous, but its derivative is not. Fortunately, the 
discontinuity in the derivative is not arbitrary, and can be calculated. To generate 
the primer across a sphere of influence: 



i oil T 


s ni 


( 2 . 11 ) 


The + indicates a value immediately after crossing the boundary, and the - indicates 
the value immediately before crossing the boundary. S is a 3X3 matrix 4 , dependent 
on the position and velocity of the sphere of influence with respect to the sun, and 
the velocity of the spacecraft with respect to the sphere. The S matrix is also 
dependent on whether the spacecraft is entering or leaving the sphere. With this in 
mind, the overall state transition matrix can be formed: 


\j/(t f ,t a ) = <>(t t , t 2 )W(t 2 )<})(t2^ t l)W(t 1 )<>(t 1 ,t ( )) (2.12) 
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where ti is the time when the spacecraft leaves the first SOI and t2 is the time when 
the spacecraft enters the second SOI. (Ktf^) and <t>(ti,to) are the state transition 
matrices for the travel in the second and the first SOI's, respectively, and <|)(t 2 ,ti) is 
the state transition matrix for travel between these two SOI's. If the trajectory 
consists of a launch, followed by an unpowered swingby and then an orbit insertion, 
the overall state transition matrix must be adjusted. This means a total of four 
extra matrices need to be multiplied into \|/(tf r t(j); one state transition matrix to get 
from the previous SOI to the swingby SOI, one W matrix to get in the SOI, another 
state transition matrix to get through the SOI, and one more W matrix to get out of 
the SOL 

With the overall state transition matrix and the initial conditions on the 
primer and its derivative known, the primer can be generated by equation 2.7. 

When this is complete, the primer history can be analyzed to find if there are 
indications that a DSM would be beneficial. By propagating the primer with the 
state transition matrices, the primer is guaranteed to be continuous. Its derivative 
will also be continuous everywhere except at the boundaries of the SOI's. This 
discontinuity is only a result of the mathematics involved in changing coordinate 
systems and the center of attraction governing the motion, and is not an indication 
of a non-optimum trajectory. The primer has been defined as a unit vector along 
the AV vectors, so the first and second necessary conditions stated earlier are 
satisfied. Since a DSM has not yet been added, the fourth condition is also not 
satisfied. Only the third condition is left to be satisfied. If the maximum magnitude 
of the primer along its history is greater than one, the trajectory is not optimal, and a 
DSM will be beneficial. To first order, the best time for the DSM is at the time when 
the primer magnitude is greatest 3 . Some representative optimal and non-optimal 
primer magnitude histories are show in Figures 1.2 - 1.5. 
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Approximating the Location of a DSM 3 

The state transition matrix 0(t-j , Iq) defines how the state of a spacecraft at 
time tj will be affected by changes in the state at time ty: 


<Kh't () ) = 


3f(t!) 


mi) 

dfr 


a?(t t ) 




5v(t-[) 
Dvr 


<hi jjhg 

< j > 2 i ! §22 


(2.13) 


To find the location of a DSM, consider perturbations about the trajectory to 
which a DSM is being added. The nominal and perturbed trajectories are shown in 
figure 2.1. Take the point on this nominal trajectory where the primer has its 
maximum magnitude (call this time tjn)/ and perturb the location by some amount, 
dr(t m ). From the definition of the state transition matrix: 


df(tm) 

3v(t m ) 


= <t>(t m ,t()) 


ar 0 

av 0 


(2.14) 


a?(t m )' 

3v(t m ) 


= <Ktm/tf) 


a?(t f ) 

9v(t f ) 


(2.15) 


The endpoints on the trajectory are fixed, and so 5r(tf) = 5r(ty) = 0 . Also note 
that 3r(t m )_ = 9r(t m ) + for continuity. 

From equations 2.14 and 2.15: 


9v(t m ) + -c)v(t m L = <h 2 (W t f)dv(tf)-<t> 22 (tm'h))dv(t () ) (2.16) 


Dv(t f ) = <}) 1 2 1 (t nv tf)a?(t rn ) + (2.17) 

dv(ty) = 012 _1 (tnri,to)^(tm)- ^ 2 ‘ 18 ^ 


12 


9 v(t m ) + -dv(t m )_ = [<}>22( t nv t f) < J ) 12 ^m/to^rftm) 


or 


where 


dr(t m ) = A _1 [9v(t m ) + -9v(t m )_] 


(2.19) 


( 2 . 20 ) 


A = [<))22 dm/ l f Hi 2 1 ( t m/ t f) _< l ) 22( t nv t ()) ( l ) 12 ^nvh))] ^ 2 - 21 ^ 

The thrust vector of the DSM must be aligned with the primer vector on the new 
trajectory. For small perturbations, the direction of the primer on the new trajectory 
can be taken to be the same as for the nominal trajectory. This yields the equation: 


9r(t m ) = dA 


-l Wtm) 


|W‘m> 


( 2 . 22 ) 


where d is the magnitude of the DSM AV. Jezewski and Rozendaal 3 expressed the 
total AV on the trajectory as a second order function of d. Taking the derivative of 
this expression with respect to d, and setting the resulting expression equal to 0, they 
arrived at an expression for d to be used in finding 9r(t m ). The expression can be 
added to the radius on the nominal trajectory to arrive at the position of the DSM: 


?DSM = fNOM( t m) + ^( t m) (2.23) 

A non-linear optimizer is then invoked with the initial estimates for the 
three coordinates of the DSM and tm as four independent variables, and total AV as 
a cost function. The original trajectory now consists of two legs: one from to to tm/ 
and the second from tm to If. The primer vector can then be found along each of 
these new legs, and the optimality conditions can be checked again. 


13 


Finding DSM's on Trajectories with Swingbv's 

The problem previously mentioned in incorporating a DSM either before or 
after a swingby has been solved by using the method of the previous section. With 
an initial estimate for the location and time of the DSM, the swingby parameters can 
be calculated, and the total AV for the trajectory can then be found. Two options are 
available to calculate the primer along a trajectory which uses a powered swingby. 
The first is to treat both the powered and unpowered swingby scenarios similarly. 
This would involve propagating the primer vector all the way through the SOI of 
the swingby planet. For a powered swingby where the impulse takes place at 
periapsis of the inbound hyperbolic trajectory, one state transition matrix for the 
inbound hyperbolic orbit, and a second state transition matrix for the outbound 
hyperbolic orbit are needed. As mentioned previously, an unpowered swingby 
involves four additional matrix multiplications in the calculation of \j/(tf,tQ)- A 
powered swingby would therefore require five multiplications. This method 
guarantees continuity of the primer and its derivative at the swingby, but it does not 
guarantee that the primer will be a unit vector aligned with the thrust at the point 
when the impulse for the swingby is applied. 

The second option is to analyze the trajectory in two separate parts. The first 
part is from the beginning of the leg to the swingby, and the second is from the 
swingby to the end of the leg. This method can only be used on powered swingby's, 
since an unpowered swingby does not provide a thrust vector which determines the 
boundary conditions on the primer vector. This method guarantees that the primer 
will be a unit vector aligned with the thrust at the swingby (since it will be defined 
as such), but it does not guarantee continuity of the derivative of the primer at the 
swingby. After a comparison of these two methods, which will be discussed later, 
the second option was selected. 
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3. Software Developed 


Integration into IPREP 

The methods described earlier to find the time and location of a DSM along 
an interplanetary trajectory were implemented in a computer program. The 
program is written in FORTRAN, and was originally written to run independently 
of any other software. It has recently been integrated into IPREP. Presently, it takes a 
trajectory from IPREP, finds the overall state transition matrix for each leg of the 
trajectory, and computes the primer vector along the entire trajectory. Then, the 
primer vector history is examined to find if and when the magnitude exceeded 
unity. If it has, a routine is called to generate the initial estimate for the location of a 
DSM, and this estimate along with the time of the DSM are sent to an optimizer, 
where total AV is the objective function. When the optimization is complete, the 
new trajectory is compared to the original to verify that the cost function has been 
reduced by the DSM. 

A considerable amount of work had to be performed in order to calculate the 
overall state transition matrix with the information in IPREP. As noted previously, 
with the technique used in IPREP, the magnitudes of the required AV's are found, 
but the direction of the burns is never calculated. This information is necessary to 
provide the boundary conditions on the primer vector, and therefore must be 
extracted. 

For a launch, the direction of the AV is the same as the direction of the 
velocity vector at periapsis, while for an orbit insertion, the direction of the AV is 
the opposite as the direction of the velocity vector at periapsis. In order to 
determine the direction of the AV, the part of the trajectory which takes place within 
the SOI must be taken into account. The turn angle which is incurred by the 
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velocity vector while within the SOI is one half of the turn angle of a swingby with 
the same Voo vector and periapsis radius, and can therefore be found with equation 
2.2. So the direction of the periapsis velocity can be found by rotating the Voo vector 
in the plane of the hyperbolic orbit by this turn angle. The plane of the hyperbolic 
orbit must be defined. The inclination of the orbit is taken to be equal to the 
declination of the Voo vector. As before, any plane change necessary to comply with 
input restrictions on the parking orbit will be treated as a separate burn. Then, the 
line of nodes is taken to be perpendicular to the Voo vector. The line of nodes and 
the inclination completely define the plane of the orbit, and there is now enough 
information to determine the direction of the burn. The direction of a swingby AV 
can be found in a similar manner for both the optimized and the simple swingby 
methods. 

Optimizer 

Two optimizers were implemented in the software to determine which 
optimizer was the better in terms of speed and finding a lower total AV. The first 
was a simple, first order optimizer which was written explicitly for this application. 
There are four independent variables; the three coordinates for the location of the 
DSM, and the time of the DSM. The initial estimates for these variables, which can 
be obtained as described earlier, are input to the optimizer. The original two burn 
trajectory is now broken into two separate trajectories. The Lambert problem is 
solved from the location of the first burn to the location of the DSM, and a second 
Lambert problem is solved from the location of the DSM to the location of the 
second burn. With the required transfer trajectories now known, the three AV's can 
be found, and the sum of them is saved. The optimizer then increments each of the 
three coordinates of the location of the DSM by a pre-defined step size. For the first 
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increments, this step size is taken to be one tenth of the magnitude of the position 
vector on the original two burn trajectory at the time of the DSM. The X, Y, and Z 
components of the location are incremented in both the positive and negative 
directions, and for each case, the sum of the AV's is found. The new location of the 
DSM is then assigned to be whichever of these six steps produced the greatest 
decrease in the sum of the AV's. If none of the six steps produced a reduction in AV, 
then the step size is halved for the next iteration. 

After the six steps for the location of the DSM have been tried, and the one 
with the best result is taken, the time of the DSM is incremented. The initial step 
size for the time increment is one tenth of the time of flight of the original 
trajectory. As with the other variables, the time step is applied both in the positive 
and negative directions, and the one which results in the greater decrease in AV is 
assigned as the new time for the DSM. If neither results in a decrease in AV, then 
the time step is halved, and the optimization continues with the location variables. 

Upper and lower bounds are applied on the time variable, so that the time 
step will not place the DSM outside of the time of flight for the original trajectory. 
No boundaries are put on the location variables. If one of the steps puts the location 
of the DSM in an unfavorable location, the AV will increase, and the optimizer will 
not take that step. The optimization is considered complete when the location step 
size is 0.01 of its original size, and the time step is less than 0.1 days. 

The second optimizer which was considered was a commercially available 
code known as ADS (Automated Design Synthesis) 11 . ADS offers a variety of 
optimization techniques including a conjugate gradient method and a variable 
metric technique. ADS allows the user to either supply gradient information, or 
ADS can calculate gradients internally. For this application, all the gradients were 
calculated internally by ADS, which uses a finite difference method. It was hoped 
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that ADS would converge on a minimum faster than the other, less sophisticated 
optimizer, but this was not found to be the case. The amount of time ADS required 
to complete the optimization was dependent on the method of optimization used in 
ADS. Most of the methods took significantly longer than the first optimizer. The 
method in ADS which converged as fast as the other optimizer found a solution 
which did not have as low a AV as the first optimizer. For this reason, the first 
optimizer is the one that was selected for the software. 


18 


l»il I ii l min' i Mur mu it limn imuwiimminiaiiii 



4. Results 


Cycler Trajectory Analysis 

Once the software had been written, a verification that it was working 
properly was desired. This verification was done by use of a study which had been 
done previously 12 . This study incorporated DSM's on an interplanetary cycler 
trajectory with the use of the DSM routine currently in IPREP. The cycler in this 
study is a spacecraft which is put into an orbit which continually cycles past both 
Earth and Mars. At each encounter, a swingby is performed to send the cycler back 
to the other planet. The cycler which was analyzed had eight Earth swingby's and 
seven Mars swingby's between the years 1996 and 2011. 

The results from this previous study were reproduced using the DSM routine 
(BPLANE) which is currently in IPREP. The same mission was then analyzed with 
the new DSM routine (PRIMER). When this trajectory was analyzed using the first 
option described earlier for trajectories with swingby's, a problem occurred. The 
additional matrix multiplications involved in each of the swingby's led to round off 
error in the overall state transition matrix. The propagation of the primer vector is 
very sensitive to the initial conditions, and the error in the overall state transition 
matrix caused the primer to fail to meet the final boundary conditions. For this 
reason, the second option (splitting the trajectory at each swingby) was seen to be the 
better. 

Table 4.1 shows the AV requirements for the nominal mission (where no 
DSM's are used). The AV's on legs 6 through 10 are significantly greater than those 
on the other legs, and so the analysis performed in Reference 11 involved looking 
for DSM's only on legs 6, 8, and 10. The last two columns in Table 4.1 show the AV 
requirements for the trajectories which were found with the BPLANE routine and 
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the PRIMER routine. Each was able to locate DSM's which significantly reduce the 
AV requirements of the mission. The fact that the new PRIMER routine found a 
trajectory with a slightly lower total AV is not as significant as that it did so in 
approximately one third the time. 

A closer look at leg 6 proves revealing. The DSM which was found by the 
BPLANE routine reduces the sum of the AV's for the two swingby's surrounding 
that leg to 0.562 km/sec. The total AV for the same leg on the PRIMER trajectory is 
0.624 km/sec. On this particular leg, the BPLANE routine found a more efficient 
DSM. Figure 4.1 reveals that the initial estimate for the time of the DSM made by 
the PRIMER routine was July 29, 2002. Figure 4.2 then shows that the optimizer 
moved this date back to July 23, 2002. However, the optimizer was unable to find 
the more efficient DSM on May 25, 2002 which was found by the BPLANE routine. 
The conclusion can be drawn that each routine is susceptible to local minima in the 
optimization process. Figures 4.3 to 4.6 show the primer magnitude histories for 
legs 8 and 10. 

A second analysis was performed on the same mission. Using the same 
encounter dates for each of the swingby's as the previous analysis, DSM's were 
sought on each leg. The results for this study are in Table 4.2. Surprisingly, the 
BPLANE routine finds a mission with a higher total AV than the BPLANE routine 
with three DSM's. Most of this increase in AV occurs at the Mars swingby before leg 
2. The reason that the BPLANE routine places a DSM on this leg is that the 
optimization in the BPLANE routine cannot account for the AV of a swingby which 
follows a DSM. The BPLANE routine is invoked in IPREP as soon as there is a leg 
which calls for a DSM. If this leg ends with a swingby, the AV cannot be found since 
the next leg has not yet been found to determine the outbound Voo vector for the 
swingby. This problem does not occur if the leg ends with an orbit insertion. The 
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PRIMER routine lets IPREP calculate the entire trajectory before it is invoked. Then, 
with all of the AV's known, the primer analysis can be done, and the DSM's are 
inserted last. The PRIMER routine is able to find a total of 7 DSM's on the cycler 
trajectory and reduces the total AV to 2.212 km /sec. 

Launch Window Analysis 

Another analysis involved the determination of the initial mass of a 
spacecraft which was to perform a direct (no swingby's), round trip mission to Mars. 
The initial mass which was needed to perform the mission was plotted against the 
launch date of the mission in Figure 4.7.1. Launches were taken every 5 days over a 
period of 500 days, starting from January 1, 2010. Outbound and inbound flight 
times were varied for each launch date in order to find the combination of 
encounter dates which produce a minimum initial mass. The Mars stay time was 
fixed at 60 days for all trajectories. The parameters of the mission are summarized 
in Table 4.3 and were the same for cases with or without DSM's. Parking orbit 
inclinations were not specified, and so IPREP always selects the inclination to equal 
the declination of the inbound or outbound Voo vector for orbit insertion or launch, 
respectively. 

Figure 4.7.1 shows noticeable improvements by using DSM's over the second 
half of the time period examined. Smaller improvements were also made for 
launch dates in January and February of 2010. Of the 101 launch dates examined, 98 
were improved by the use of DSM's. Figure 4.7.2 shows the percentage of mass on 
the nominal mission which is saved by the use of DSM's. The discontinuity on this 
plot is caused by the three dates where no DSM could be found, since the zero 
percent savings cannot be plotted on the logarithmic scale. The purpose of this 
study was not to do a detailed examination of when launch opportunities occur, but 
rather to show that DSM's can expand the launch windows. Figure 4.7.1 shows that 


21 



by using DSM's, the minimum initial mass pertaining to a specific launch date may 
be reduced below a specified maximum allowable initial mass, thus including that 
date in a launch window. 

Abort Mission Analysis 

Two round trip Mars missions were selected to serve as nominal cases for an 
abort mission analysis. These missions were selected because they have a relatively 
low total AV and initial mass, and could be considered practical candidates for a 
manned Mars mission. Again, this study was done to show the usefulness of the 
DSM capability in IPREP, and not as an in-depth mission analysis. 

The launch date of the first mission is in the year 2020, and the mission 
parameters are outlined in Table 4.4. The nominal mission (where there is no 
abort) has a two year trip time, including 60 days stay time at Mars. The abort 
mission assumes that some problem occurred after launch, and the spacecraft is to 
return to Earth at least 50 days earlier than the nominal mission. To accomplish 
this, the orbit insert at Mars is replaced with a swingby. The Mars encounter date for 
the abort mission is the same as for the nominal mission, while the Earth return 
date is allowed to vary up to the maximum trip time to find the minimum initial 
mass needed to complete the mission. The abort mission was run with and without 
searching for a DSM on the return leg. A DSM was not used on the first leg, since it 
is assumed that the problem occurred after launch, and the vehicle has achieved the 
desired trajectory which will rendezvous with Mars. The second mission has a 
launch date in the year 2022. Its parameters are also outlined in Table 4.4, and the 
abort scenario is the same as for the first case. 

The desired result is for the abort trajectory to require less initial mass than 
the nominal case. When this is true, the abort trajectory is available at no additional 
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cost. If the abort trajectory requires a slightly greater initial mass, then it may still be 
practical to include this extra propellant mass on the nominal mission in order to 
have the abort trajectory available should a problem occur. The results of this study 
are shown in Tables 4.5.1 and 4.5.2. Unfortunately, a practical abort trajectory was 
not available on either of these two missions. However, DSM's were able to 
improve the abort trajectories in each case. Although the improvement was not 
great enough to make these abort trajectories feasible, the indication is that for 
future abort mission analyses, it is worthwhile to search for DSM's. Other abort 
studies have tried using Venus swingby's to open up abort trajectories. Having the 
DSM capability in IPREP creates more possibilities for abort analysis such as using 
Venus swingby's in combination with DSM's. 

Figures 4.8 and 4.9 show the primer magnitude history for the 2020 abort 
mission with and without the DSM, respectively. Figures 4.10 and 4.11 are similar 
figures for the 2022 abort missions. The trajectories for the nominal and abort 
missions for both the 2020 and the 2022 missions are shown in figures 4.12 to 4.17. 
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5. Conclusions 


Primer vector theory is useful in determining if a trajectory will be improved 
by the use of a DSM. Once the primer magnitude history has determined that a 
DSM should be incorporated, primer vector theory also provides a method to 
determine initial estimates for when and where to place the DSM. The problem of 
being unable to calculate swingby parameters when a swingby is either preceded or 
followed by a DSM is avoided by the use of primer vector theory. This allows any 
trajectory to be subject to a primer analysis. 

The software developed was shown to work effectively and efficiently in 
improving interplanetary trajectories with deep space maneuvers. 

The cycler trajectory analysis verified that the software developed will find 
DSM's which reduce the AV required on a trajectory. The results compared 
favorably to the BPLANE routine in IPREP which also finds DSM's. The greatest 
advantage of the new software is the reduction in time required to find DSM's. 
Although each routine found trajectories with comparable AV's, the dates of the 
DSM's vary considerably, indicating that each code is susceptible to local minima in 
the optimization of the DSM's. However, the new software is more robust, in that it 
will never increase the total AV on a trajectory, while the BPLANE code may do so. 

The launch window analysis showed another use for DSM's. DSM's were 
shown reduce a vehicle's initial mass required to complete a mission. If a 
maximum initial mass is defined, it is apparent that DSM's may be able to reduce 
the initial mass associated with some launch dates from above this maximum to 
below it, thus expanding the launch windows for that mission. 

The final study showed that DSM's are also useful in looking for possible 
abort trajectories. Although no practical abort trajectories were found for the two 
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mission examined, the data indicated that some abort trajectories will become 
available with the use of DSM's. 
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6. Recommendations for Future Work 

The software developed can be used to perform more detailed analyses of the 
type which were done for this paper. Launch window analyses for direct Mars 
missions as well as missions which use Venus swingby's on the inbound or 
outbound leg can be performed. More in-depth abort analyses can be performed. 
'Grand Tour' missions which pass by many planets on the way out of the solar 
system can be done with DSM's on some or all of the legs. 

As for more development in the software, the most important suggestion 
would be to have the capability to look for more than one DSM on one leg of a 
trajectory. The software to propagate the primer vector along a trajectory after a 
DSM has been inserted has already been written. With the primer magnitude 
history known, all that remains to be done is to generate an initial estimate for the 
location of the second DSM, and then to invoke an optimizer. The optimizer would 
have to be altered for the second DSM. The time and three position coordinates of 
both the first and second DSM's could be used as independent variables in the 
optimization, and the optimizer currently in use was written specifically for the 
optimization of one DSM with four independent variables. Figures 4.2, 4.4, and 4.6 
indicate that a second DSM could be used on those legs to further reduce the AV. 

Non-linear methods could be applied to the problem when primer vector 
theory fails due to the limits imposed by the linear approximations. Non-linear 
methods could improve the initial estimate for the location of the DSM when that 
location is distant from the nominal trajectory, and make the software more robust. 


26 



7. References 


1. Fimple, W.R., "Optimal Midcourse Plane Change for Ballistic Interplanetary 
Trajectories, AIAA Journal, Vol. 1, No. 2, Jan. 1963, pp. 430-434. 

2. Lawden, D.F., Optimal Trajectories for Space Navigation, Butterworths, 
London, 1963, pp. 3-25, 54-66. 

3. Jezewski, D.J., and Rozendaal, H.L., "An Efficient Method for Calculating 
Optimal Free-Space N-Impulse Trajectories," AIAA Journal, Vol. 6, No. 11, 
Nov. 1968, pp. 2160-2165. 

4. Glandorf, David R., "Primer Vector Theory for Matched-Conic Trajectories," 
AIAA Journal, Vol. 8, No. 1, January 1970, pp. 155-156. 

5. Lion, P.M., and Handelsman, M., "Primer Vector on Fixed-Time Impulsive 
Trajectories, AIAA Journal, Vol. 6, No. 1, Jan. 1968, pp. 127-132 

6. Friedlander, Alan L., "MULIMP - Multi-Impulse Trajectory and Mass 
Optimization Program", Report No. SAI 1-120-383-T4, Science Applications, 
Inc., for Jet Propulsion Laboratory, April 18, 1975. 

7. Hong, P.E., Kent, P.D., Olson, D.W., and Vallado, C.A., "Interplanetary 
Program to Optimize Simulated Trajectories (IPOST) User's Guide", Martin 
Marietta Astronautics, Space Launch Systems Company, Denver, Colorado, 
NASA CR-1 89653, Volume I, October, 1992. 

8. Hong, P.E., Kent, P.D., Olson, D.W., and Vallado, C.A., "Interplanetary 
Program to Optimize Simulated Trajectories (IPOST) Analytic Manual", 
Martin Marietta Astronautics, Space Launch Systems Company, Denver, 
Colorado, NASA CR-189653, Volume II, October, 1992. 

9. Striepe, Scott Allen, "Effect of Interplanetary Trajectory Options on Entry 
Velocities at the Earth and Mars Atmospheric Interface," Thesis, May 1991, 
The University of Texas at Austin. 

10. Battin, Richard FI., An Introduction to the Mathematics and Methods of 
Astrodynamics, Second Printing, AIAA, Inc. New York, NY, 1987, pp. 450-463. 

11. Vanderplaats, G.N., "ADS - A Fortran Program for Automated Design 
Synthesis - Version 1.10", NASA Contractor Report 177985, Sep. 1985. 

12. Fruth, Gregory, "Establishing an Aldrin Cycler Between Earth and Mars Using 
Nuclear Electric Propulsion," Thesis, February, 1993, The George Washington 
University, Washington, D.C. 


27 


13. Baker, Jerome M., "Orbit Transfer and Rendezvous Maneuvers between 
Inclined Circular Orbits", Journal of Guidance, Control and Dynamics, Vol. 3, 
No. 8, Aug. 1966, pp. 1216-1220. 

14. Bate, Roger, Mueller, Donald, and White, Jerry, Fundamentals of 
Astrodynamics, Dover Publications, 1971. 

15. Bryson, A.E. Jr., Denham, W.F., and Dreyfus, S.E., "Optimal Programming 
Problems with Inequality Constraints, I: Necessary Conditions for Extremal 
Solutions", AIAA Journal, Nov. 1963. 

16. Bryson, Arthur E., and Ho, Yu-Chi, Applied Optimal Control , Hemisphere 
Publishing Corp., New York, NY, 1975. 

17. Eckel, Karl, "Numerical Solutions of Non-Coaxial Optimum Transfer 
Problems, The Journal of the Astronautical Sciences, Vol. X, No. 3, pp. 82-92, 
Fall 1963. 

18. Edelbaum, T.N. "How Many Impulses?", AIAA Paper No. 66-7, AIAA 3rd 
Aerospace Sciences Meeting, New York, Jan. 24-26, 1966. 

19. Edelbaum, T.N., "Minimum Impulse Transfers in the Near Vicinity of a 
Circular Orbit", The Journal of Astronautical Sciences, Vol. XIV, No. 2, pp. 
66-73, Mar.-Apr. 1967. 

20. Gerald, Curtis F., and Wheatley, Patrick O., Applied Numerical Analysis, 
Fourth Edition, Addison-Wesley Publishing Company, New York, 1989. 

21. Glandorf, David R. "Lagrange Multipliers and the State Transition Matrix for 
Coasting Arcs" AIAA Journal, Technical Notes, Volume 7, Number 2, Feb. 
1969, pp. 363-365. 

22. Hiller, H. "Optimum Impulsive Transfers Between Elliptic and Non- 
Coplanar Circular Orbits", Planet, Space, Science, Vol. 13, pp. 1233-1247, 
Programming Press, Ltd., 1965. 

23. Hiller, H. "Optimum Impulsive Transfers Between Non-Coplanar Elliptic 
Orbits Having Colinear Major Axes", Planet, Space, Science, Vol. 14, pp. 773- 
789, Programming Press, Ltd., 1965. 

24. Hiller, H. "Optimum Transfers Between Non-Coplanar Circular Orbits", 
Planet, Space, Science, Vol. 13, pp. 147-161, Programming Press, Ltd., 1965. 


28 


'limn IN" in Itlrfil h il Mil li n i ii 1 111111.11 lii< jin jnnr* iwim m « iHm il kiIimI'IH nil 1111 rvttmin Mli mill 11 mil mm 1 1 1 mi ft I Mill wiiMiiiiiiiMMiwiiii lai'riHMHiMMHmmiiNini 


25. Jezewski,. Donald J., "Primer Vector Theory and Applications", NASA TR R- 
454, Nov. 1975. 

26. Kirk, Donald E., Optimal Control Theory, An Introduction, Prentice Hall, 
New Jersey, 1970. 

27. Marec, J.P. Optimal Space Trajectories, Elsevier Scientific Publishing 
Company, New York, 1979. 

28. Meirovitch, Leonard, Methods of Analytical Dynamics, McGraw-Hill, 1970. 

29. Miele, A., Wang, T., and Basupar, V.K., "Primal and Dual Formulations of 
Sequential Gradient-Restoration Algorithms for Trajectory Optimization 
Problems", Acta Astonauiica, Vol. 13, No. 8, pp. 491-505, 1986. 

30. Prussing, John E., and Chiu, Jeng-Hua, "Optimal Multi-Impulse Time-Fixed 
Rendezvous Between Circular Orbits", Journal of Guidance, Control and 
Dynamics, Vol. 9, No. 1, Jan. -Feb. 1986, pp. 17-22. 


29 


8. Tables 


Table 4.1: Comparison of BPLANE Code and PRIMER Code for Cycler Trajectory 

from Reference 12 



NOMINAL 


LEG 

NUMBER 


SWINGBY 

PLANET 


DATE 


AV 

(KM/SEC) 


BPLANE 


AV 

(KM/SEC) 


0.073 


0.081 


0.067 


0.064 


0.098 


0.000 




PRIMER 


AV 

(KM/SEC) 


0.073 


0.081 


0.067 


0.064 


0.098 


0.004 





EARTH 


MARS 


EARTH 


MARS 


EARTH 


MARS 


DSM 


EARTH 


MARS 


DSM 


EARTH 


MARS 


DSM 


EARTH 


MARS 


EARTH 


MARS 


EARTH 


19 NOV 1996 


1 MAY 1997 


1 JAN 1999 


28 MAY 1999 


8 FEB 2001 


6 JUL 2001 


25 MAY 2002 


23 JUL 2002 


16 APR 2003 


12 SEP 2003 


9 AUG 2004 


15 JUL 2004 


7 JUL 2005 


13 DEC 2005 


17 OCT 2006 


27 JUL 2006 


6 SEP 2007 


16 FEB 2008 


10 OCT 2009 


28 MAR 2010 


13 NOV 2011 


0.067 


0.064 


0.098 


1.536 


2.467 


2.247 


3.699 


1.062 


0.261 


0.105 


0.074 


0.125 


0.110 



0.100 


0.000 




0.116 


0.000 




0.126 


0.105 


0.074 


0.125 


0.110 


0.488 


0.132 


0.001 



0.619 


0.042 


0.002 



0.419 


0.023 


0.105 


0.074 


0.125 


0.110 


TOTAL AV 
(KM/SEC) 

12.069 

2.565 

2.528 

RUN TIME 
(SEC) 

3.7 

53.2 

19.1 
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Table 4.2: Comparison Between BPLANE Code and PRIMER Code for Cycler 
Trajectory When a DSM is Allowed on Each Leg 



BPLANE 

PRIMER 

LEG 

NUMBER 

SWINGBY 

PLANET 

DATE 

AV 

(KM/SEC) 

AV 

(KM/SEC) 

1 

EARTH 

19 NOV 1996 

0.000 

0.073 

DSM 

20 DEC 1996 

0.072 





2 

MARS 

1 MAY 1997 

1.501 

0.003 

DSM 




27 MAY 1998 


0.041 

3 

EARTH 

1 JAN 1999 

0.645 

0.027 

DSM 


• , ' : • || pi 1 1 





4 

MARS 

28 MAY 1999 

0.000 

0.002 

DSM 

26 JUL 2000 

0.030 


23 MAY 2000 


0.060 

5 

EARTH 

8 FEB 2001 

0.000 

0.070 

DSM 

8 MAR 2001 

0.103 





6 

MARS 

6 JUL 2001 

0.000 

0.004 

DSM 

25 MAY 2002 

0.457 


23 JUL 2002 


0.488 

7 

EARTH 

16 APR 2003 


0.132 

DSM 

18 MAY 2003 






8 

MARS 

12 SEP 2003 

0.000 

0.001 

DSM 

9 AUG 2004 

0.577 


15 JUL 2004 


0.619 


continued 
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Table 4.2 Concluded 



BPLANE 

PRIMER 

LEG 

NUMBER 

SWINGBY 

PLANET 

DATE 

AV 

(KM/SEC) 

AV 

(KM/SEC) 

9 

EARTH 

7 JUL 2005 

0.000 

0.042 


DSM 

5 AUG 2005 

0.225 





10 

MARS 

13 DEC 2005 

0.000 

0.002 

DSM 

17 OCT 2006 

0.370 


27 JUL 2006 


0.419 

11 

EARTH 

6 SEP 2007 

0.000 

0.023 

DSM 

30 SEP 2007 

0.188 





12 

MARS 

16 FEB 2008 

0.000 

0.001 

DSM 

25 MAR 2009 

0.071 


15 FEB 2009 


0.075 

13 

EARTH 

10 OCT 2009 

0.000 

0.009 

DSM 

7 NOV 2009 

0.012 





14 

MARS 

28 MAR 2010 

6.000 

0.003 

DSM 

14 MAY 2011 

0.079 


11 MAR 2011 


0.069 

EARTH 

13 NOV 2011 

0.039 

0.048 


TOTAL AV 
(KM/SEC) 

4.533 

2.212 

RUN TIME 
(SEC) 

332.2 

54.6 
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Table 4.3: Vehicle and Mission Parameters for Round Trip Mars Mission Launch 

Window Analysis 


Vehicle 

Mars drop off mass, kg 76000 

Earth return mass, kg 61000 

Engine Specific Impulse (Isp), sec 480 


Mission 

Earth departure parking orbit 

Mars parking orbit 

Earth return parking orbit 

Outbound time of flight 

Mars stay time 

Inbound time of flight 

Maximum total time of flight 


..500 km altitude, circular 

..500 km periapsis altitude, 1 sol 

..500 km periapsis altitude, 1 sol 

.108 to 418 days 

.60 days 

.73 to 383 days 

730 days 
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Table 4.4: Vehicle and Mission Parameters for Nominal Mars Mission Abort 

Analysis 


Vehicle 

Mars drop off mass, kg 76000 

Earth return mass, kg ..61000 

Engine Specific Impulse (Isp), sec 475 


Mission 

Earth departure parking orbit 500 km altitude, circular 

Mars parking orbit 500 km periapsis altitude, 1 sol 

Earth return parking orbit 500 km periapsis altitude, 1 sol 

Dates 



2020 Mission 

2022 Mission 

Launch Date 


11 Oct. 2022 

Outbound time of 
flight 

278 days 

288 days 

Mars stay time 

60 days 

60 days 

Inbound time of 
flight 

383 days 

383 days 

Total time of flight 

721 days 

731 days 
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Table 4.5.1 Comparison Between Nominal and Abort Trajectories 

for 2020 Mars Mission 



Nominal 

Abort 

Without 

DSM 

Abort 
With DSM 

Initial Mass 
(million kg) 

1.6 

9.9 

4.1 

Total AV 
(km /sec) 

14.5 

23.6 

19.4 

Time of flight 
(days) 

721 

651 

501 


Table 4.5.2 Comparison Between Nominal and Abort Trajectories 

for 2022 Mars Mission 



Nominal 

Abort 

Without 

DSM 

Abort 
With DSM 

Initial Mass 
(million kg) 

2.7 

4.5 

3.8 

Total AV 
(km /sec) 

15.0 

19.8 

18.9 

Time of flight 
(days) 

731 

681 

511 
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Figfure 2.1: Perturbations from the Nominal Trajectory 
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1 . Earth Departures 1 Sep 2020 

2. Mars Arrival: 16Jun 2021 

3. Mars Departure: 15 Aug 2021 

4. Earth Arrival: 2 Sep 2022 


f igure 4.12: 2020 Round Trip Mars Mission 




1. Earth Departure: 11 Sep 2020 

2. Mars Swingby: 1 6 Jun 2020 

3. Earth Arrival: 24 Jun 2022 


Figure 4.13: 2020 Abort Mars Mission Without DSM 
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1 . Earth Departure^ 1 Sep 2020 

2. Mars Arrival: 16 Jun 2021 

3. DSM : 1 Nov 2021 

4. Earth Arrival: 25 Jan 2022 


Figure 4.14: 2020 Abort Mars Mission With DSM 



1. Earth Departure: 1 1 Oct 2022 

2. Mars Arrival: 26 Jul 2023 

3. Mars Departure: 24 Sep 2023 

4. Earth Arrival: 11 Oct 2024 


Figure 4. 1 S: 2022 Round Trip Mars Mission 
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1. Earth Departure: 11 Oct 2022 

2. Mars Swingby: 26 Jul 2023 

3. Earth Arrival: 22 Aug 2024 




t 




Figure 4.16: 2022 Abori Mars Mission Without DSM 



1. Earth Departure: 11 Oct 2022 

2. Mars Arrival: 26 Jul 2023 

3. DSM 7 Dec 2023 

4. Earth Arrival: 5 Mar 2024 


Figure 4.17: 2022 Abort Mars Mission With DSM 
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10. Appendix 


Derivation of the Necessary Conditions on the Primer Vector 

The derivation presented here is a summary of the original derivation done 
by D.F. Lawden^. The necessary conditions on the primer vector provide the basis 
for this paper, and an understanding of the derivation will give the reader a 
stronger grasp of the subject. Indicial notation is used in this appendix. A subscript 
which appears only once in any term is known as a 'free index', and can take on any 
/ value in its range. A subscript which appears twice in a term is known as a 'dummy 
index'. A dummy index implies a summation of all terms over the range of the 
index. 

Define a cost function, J(xi,x 2 ,...,xn,tf), where x is an n-dimensional state 
vector. The problem is to minimize J, subject to: 


xj = fj(x,a,t) for i = 1 to n, (A.l) 

h k (x,a,t) = 0 for k = 1 to p < m. (A.2) 

Xj = XjQ for i = 1 to n at time = to (A.3) 

x w = x wf for w = 1 to q < n at time = tf (A.4) 


Here, fi are the equations of motion, hk are constraint functions, t is time, and a is 
an m-dimensional control vector. Let x*(t) and a*(t) be the state and control 
histories which satisfy the equations A.l to A.4 and also minimize J. Introduce a 
small parameter e, and use it to perturb these state and control vectors. That is, let 

x(t,e) = x*(t) + e(t) and a(t,e) = a*(t) + e(t). 

Differentiate equation A.l with respect to e to obtain 
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(A. 5) 


a 2 Xi = afj 9x r 1 3fj ^ a j 
3e3t 3x r de 3aj de 


Now, let 


yi<t) = 


5xj 

. 3e A = 


e=0 


and Pj(t) = 

so at e = 0, equation A. 5 becomes 


(d aj 
v 3e y 


e=0 


3fj 3fj 

yi= s; y ' + 3^i 


(A. 6) 


Similarly, by taking the derivative of equation A. 2 with respect to e: 


^lk y +^kp =0 (A. 7) 

3xj y ‘ 3aj 

Now the Lagrange expression is introduced: 

F = -Xjfj +p k h k (A.8) 

This is where the primer vector is introduced. The primer vector is a Lagrange 
multiplier. In problems of the type being considered here, the state vector is a 7- 
dimensional vector, with the first three components being the velocity of a 
spacecraft, the next three components being the position of the spacecraft, and the 
last component being the mass of the spacecraft. The components of the primer 
vector are, by definition, the first three components of X. It will be shown that the 
derivative of the primer vector is a vector whose components are the next three 
components of X. 

From the definition of F: 
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J 

to 


, . 3F OF _ 

hy\ + t~ + ^ — P i 


dx\ 




^ tf 

f f 

dt= J 

h 

y ‘o 

v V 


afi 5/j R 


+^k 


Oh 


dx T 


k„ 




y r + 


a a . n 

da j i) 


dt 
(A. 9) 


By equations A. 6 and A. 7, this expression equals 0 on an optimal trajectory. 
Integration by parts yields: 




where 



OF ^ 
3^ dl 


The parameters Xi and pk are at our disposal, and so Xi can be chosen such that 

Xj = Gj +X.JQ (A. 10) 

where Xjo are constants yet to be determined, and pk can be chosen such that 
dF OF 

— — = 0. Equation A. 10 can also be expressed as Xj = . Equation A. 9 can now be 

0aj< Ox; 

reduced to 


r OF 

o = yif^if-yi(Aio+ J ^-Mt 


(A.ll) 


‘0 


where r = p+l,p+2,...,m and Xjf = Xj(lf). This result will be used later. 

Define a family of variations e a , where o runs from 1 to N, and N = q + n + 1. 
Each member of this family satisfies equations A.l to A.4. Now J can be expressed as 
a function of these variations: 


J = J(e 1 ,e 2 ,...e N ) = ] 0 + U 

where Jo is the minimum J. U = 0 if and only if e°= 0 for all a. In all other cases, U 
must be positive, since Jo is by definition the minimum J. Theory of implicit 
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functions tells us that e ® can be determined as a continuous function of U in the 
neighborhood of U=0 if and only if the determinant 


A = 


JL 

aj 

a J 

ae 1 

3e 2 

ae N 

ax i0 

ax i0 

ax i0 

De 1 

ae 2 

ae N 

ax W f 

ax wf 

ax W f 

ae 1 

ae 2 

ae N 


is not zero. If the determinant is not zero, it implies that there are some variations 

V 

which produce a negative value of U. This cannot be the case, since U must be \ 
greater than or equal to zero. The conclusion is that the determinant must vanish, 
and therefore that the rows of the above matrix must be linearly dependent. This 
means that N constants yo, yj (i = 1 to n), and v w (w = 1 to q) can be found such that 


\ 


dj 3x i( ) . 0x w f n 

Y()TT + YiTr + v l“T^ = 0 

De z ae de 


(A. 12) 


for arbitrary z in the range of 1 to N. Looking at each of the three terms in this 
equation when e z = 0 (optimal trajectory): 


where 


jn _ 3J ( **, z , ^x s 
l9e z J e z =o ^ x sfv^f d£ z 

s = q+1, q+2, ... n 
x s = free states at t = tf 


3J z 

+ a^ Uf 


Uf z = 


dtl) 

dE z J e z =0 


The second and third terms become: 
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(no sum over 0 


Dx: 


de 


10 _ v 
t~ y«o 


X 


5x w f - z 
— ^ = x w fUf + y w f 


/ 

/ Inserting these three expressions into A. 12, and also adding in equation A. 11: 


( 31 

® — (Yi — ^io)yiO (^w +^wf)ywf+ Yo^ ^ *sf 

V dx sf 


y S f + 


aj . aj 
To x sf + Yo~, + v w x wf 
, ^ x sf ^t f 


‘0 r 


(A. 13) 


The ranges of the indecis in equation A. 13 are as follows: 
i = 1 to n (n = number of states) 
w = 1 to q (q = number of fixed final states) 
s = (q+1) to n 

r = (p+1) to m (p = number of constraints, m = number of controls) 

The subscripts 0 and f indicate initial and final conditions, and are not indecis to be 
summed over. 

In order for equation (A. 13) to vanish for arbitrary variations, each of the 
terms must vanish: 

X i( ) = Yi (A.14) 


Vw — ^ wf 


d) , 

Yo— = -*sf 


dx K 


sf 


aj aj . „ 

Y{) -v + Yo x sf + v w x wf “ ® 

atf ax sf 


(A.15) 
(A. 16) 

(A. 17) 
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(A. 18) 


_3F 
da r 


= 0 


Substituting (A.15) and (A.16) into (A.17), and realizing that 


^■sf*sf ^ wf*wf 


the result is 


, • 3J 


\ 

\ 

\. 

(A^19) 


These equations are now ready to be solved. There are 2n + m + p unknown \ ' 

functions which are the n state variables Xj, the m control variables aj, the n ^ 

v « 

functions X\, and the p functions pk- We have the same number of equations in 
(A.l: n equations), (A. 2: p equations), (A. 10: n equations), and (A. 18: m equations). 

We also have the n equations A. 3, q equations A. 4, (n-q) equations A.16, and 1 
equation A. 19 to determine the n constants Xio, the n constants of integration for the 
state variables, and tf. If tf is fixed for a particular problem, equation A. 19 no longer 
applies. If these necessary conditions can be satisfied for Yo * 0 , then the solution is 
said to be normal and yo can be chosen to equal 1 in equation A. 12 without loss of 
generality. 

Examine a trajectory problem defined as follows: 


cm . 

v i = 77 i + Si 
M 


M = -m 

cm 

where i = 1,2,3 for the X,Y, and Z coordinates, — is the magnitude of the engine 

M 

thrust per unit mass, lj are the direction cosines of the thrust vector, and g[ are the 
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components of gravitational acceleration. M is the mass of the spacecraft and m is 
the propellant mass flow rate. Two constraint equations can be defined: 


h 2 +^2 2 +1 3 2 “ 1 

/ m (m max -m) = 0 

The latter equation indicates that only no-thrust or maximum-thrust arcs are 
Allowed on this trajectory. The state vector is: 

j 

[vi v 2 v 3 ^ r 2 r 3 M] 

and the control vector is: 

ft h h m ] T 

The Lagrange expression is: 

F = + 8i )“ Xi+3V ' + + + ^ + 1:)2 “ 1 ) + ^2[ m ( m max " m )] 

Employing equation A. 18 and the derivative of A.10 yields: 


>> 

II 

1 

>> 

T 

(A. 20 ) 

i _ x 

(A. 21 ) 


(A. 22) 

M 2 


0 = -^X i+ 2p 1 l i 
M 

(A. 23) 

M 

(A. 24) 


for i and j = 1 to 3. The first two of these equation can be combined to form: 
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h = x. 


agi 


Equations A.20 and A.21 indicate that the primer and its derivative are continuous 

functions of time. \ 

The Weierstrass condition^ (submitted here without proof) for this problem 

can be expressed as: 


-Xy ]m > --X.jli-X. 7 Jm 
M i (M 


for all possible values of 1 and ni which satisfy the constraints. On a no-thrust arc, \ 
m= 0 , and so the condition requires that the right hand side be negative, or that: 


V 


X .7 ^ — X;lj 
' M " 

The right hand side of this expression takes its maximum value when the primer is 
aligned with the thrust vector. Thus, the condition on no- thrust arcs becomes: 



P 


where p is the magnitude of the primer vector. 

On max-thrust arcs, two cases need to be considered. First, let m* take on its 
maximum value, m m ax’ In order for the Weierstrass condition to be satisfied: 


which implies that: 


Xilj>Mr 

Xjij >p 


Since I is a unit vector, the > will never hold, and so: 

hh =P 




56 


This indicates that the thrust vector must always be aligned with the primer vector. 
For the second case, let m* take on its minimum value of 0, then it is necessary that 
the left hand side of the Weierstrass condition be positive: 


’k.y < — A.;l; => A. 7 < — p 

on. max-thrust arcs. This leads us to define the switching function k: 


K = p - \ 7 

On no-thrust arcs, it is necessary that the switching function must be less than or 
equal to zero. On max-thrust arcs, the switching function must be greater than or 
equal to zero. If impulsive thrusts are allowed in the simulation, periods of 
maximum thrust are modeled to be instantaneous, and the switching function must 
equal zero at an impulse. If an impulse is required at a point which is not at either 
end of a trajectory (a deep space maneuver), the value of the switching function 
must be negative immediately preceding and immediately following the burn. The 
derivative of the switching function has been shown 2 to be continuous, and 
therefore, at a deep space maneuver, the derivative of the switching function must 
also be zero. 

Define the cost function, J, to be the sum of the magnitudes of the AV's over 
the whole trajectory, which, by the rocket equation, can be expressed as: 


J = c In 

With this cost function defined, 


'Mg 
v M f ; 


= go Isp In 




M 


fj 


X 7 = — = — 
* M f M 
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by equation A. 16 at time = tf. Equation A. 22 shows that X 7 is constant 0^1 a no-thrust 
arc (since m=0). Therefore, on a no-thrust arc, ' 

(A25> 

At the two impulses surrounding this no-thrust leg, the switching function equals 
zero, and so p must equal 1. Also, the switching function must never be positive on 
a trajectory consisting only of impulses and no- thrust arcs. So the maximum 
magnitude of the primer vector is unity, and that only occurs at the impulses. 
Equation A.25 indicates that 



We have said that k must be zero at a deep space maneuver, therefore, the 
derivative of the magnitude of the primer vector must also be zero there. 

Four necessary conditions to be satisfied for an optimal trajectory allowing 
impulsive thrust, and whose cost function is the sum of the AV s have been 
derived: 

1) The primer vector and its first time derivative are continuous. 

2) During any impulse, the thrust vector must be aligned with the primer. 

3) The magnitude of the primer is a maximum and has a value of one during 

any impulse. 

4) The derivative of the magnitude of the primer is zero at a deep space 

maneuver. 


!a 
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